Partition Function Zeros 
for Aperiodic Systems 

Michael Baake^'^'J], Uwe Grimm^'^, and Carmelo Pisani^'"^ 

^ Department of Mathematics, University of Melbourne, 
Parkville, Victoria 3052, Australia 

^ Institut fiir Theoretische Physik, Universitat Tiibingen, 
Auf der Morgenstelle 14, 72076 Tiibingen, Germany 

^ Instituut voor Theoretische Fysica, Universiteit van Amsterdam, 
Valckenierstraat 65, 1018 XE Amsterdam, The Netherlands 

Department of Mathematics, La Trobe University, 
Bundoora, Victoria 3083, Australia 

Preprint itap-11/93-2 

The study of zeros of partition functions, initiated by Yang and Lee, provides an important 
qualitative and quantitative tool in the study of critical phenomena. This has frequently 
been used for periodic as well as hierarchical lattices. Here, we consider magnetic field 
and temperature zeros of Ising model partition functions on several aperiodic structures. 
In ID, we analyze aperiodic chains obtained from substitution rules, the most prominent 
example being the Fibonacci chain. In 2D, we focus on the tenfold symmetric triangular 
tiling which allows efficient numerical treatment by means of corner transfer matrices. 
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1 Introduction 

The study of critical phenomena by means of discrete spin systems, initiated by Lenz and 
Ising finally came of age by Onsager's spectacular solution of the 2D field- free Ising 
model [0. He could show that the 2D Ising model on the square lattice with ferromagnetic 
nearest neighbour interaction exhibits an order- disorder phase transition of second order at 
finite temperature with the magnetization as order parameter. 

The investigation of many other spin systems followed, as did the consideration of Ising- 
type models on other lattices and graphs. Some cases can still be solved exactly, compare 
P], but this is an exception. Consequently, one needs complementary methods to tackle 
questions like existence and location of critical points and estimation of critical exponents, 
especially in higher dimensions. One such technique is provided by the work of Lee and 
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Yang 0, ^, who investigated the distribution of partition function zeros in the complex 
plane, i.e., field variables or fugacity were treated as a complex parameter. Later, also the 
zeros in the complex temperature variable were studied for various systems @, 0, §], 0]- These 
not only provide information about the location of phase boundaries but also about critical 
exponents, see e.g. |§] where mainly hierarchical lattices are considered. For those, the zero 
patterns form fractal structures known as Julia sets whereas they are usually expected to lie 



on simple curves for regular lattices, at least in the isotropic case. It has been shown [10, 11 



that for anisotropic Ising models on two-dimensional regular lattices the temperature zeros 
generically fill areas in the complex plane. 

But what does "simple curves" mean? As we will show in a shortwhile, already the Ising 
model on a modulated structure in ID can result in magnetic field zeros which do lie on a 
simple curve, but occupy only a Cantor-like portion of it - even in the thermodynamic limit! 
Although the corresponding phenomenon in 2D does not typically show up (or at least not 
that we could conclude so), the investigation of partition function zeros for spin systems 
on non-periodic graphs with deflation/inflation symmetry might fill the gap between the 
relatively well studied cases of lattices and hierarchical graphs. This is the major motivation 
for the present article, where we discuss the Ising model on quasiperiodic graphs in ID and 
2D. 



2 Non-periodic ID Ising model 

Let us consider a one-dimensional chain of N Ising spins (Xj G {±1}, j = 1, 2, . . . , A^, with peri- 
odic boundary conditions (i.e., cr7v+i = o'i)- The energy of a configuration cr = ((Xi, 0-2, • • • , (^n) 
is given by 

N 

= - ^ (Tj(Tj+i + i/j (Tj . (2.1) 

i=i 

Here, we concentrate on systems where the coupling constants Jjj+i can only take two 
different values Jjj+i G {Ja, Jb} and where the magnetic field is constant, i.e. Hj = H. The 
actual distribution of the two coupling constants Ja and Jb along the chain is determined by 
an infinite word in the letters a and b which is obtained as the unique limit of certain two 



letter substitution rules, see 12 . In fact, we restrict ourselves to the Fibonacci case which 



corresponds to the substitution rule 

- : : I. 

The length of the word Wn = S{wn-i) obtained by n iterations from the initial word wq = a 
is /n+i, where the Fibonacci numbers /„ are defined by 

/O = , /l = 1 , fn+l = fn + fn-1 ■ (2.3) 
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Let us introduce the following notation 

Ka = K. = ^, h = (2.4) 

KbT kbT kbT 

and 

Za = exp{2Ka) , Zb = exp{2Kb) , w = exp{2h) . (2.5) 
The two elementary transfer matrices and now read as follows 



Ta, = (t. ■ ^„,,)-^/M ^ ■ = {w ■ z^,b)-'/' f^,b (2.6) 



and, in general, do not commute with each other. 

The recursive definition of the sequence gives rise to a recurrence formula for the transfer 
matrices 

To = Ta , Ti = Tb , Tn+i = Tn ■ T„_i , (2-7) 

where T„ denotes the transfer matrix of the chain that corresponds to the n-th iteration step. 
Hence the partition function Zn{za, Zb, w) = tr(T„) is essentially (i.e., up to an overall factor) 
a polynomial in its three variables. It is also possible to write down a recurrence relation for 
the partition function itself, see [|l^ for details. In this way, it is really easy to generate the 



partition function for very large systems exactly (e.g., by means of algebraic manipulation 
packages). Our problem at hand has thus been reduced to the task of computing the roots 
of a polynomial. 

Let us have a look at the pattern of zeros of the partition functions Zn{za, Zb, w) in the 
field variable w. If both couplings are ferromagnetic [za > ^, Zb > 1), the zeros lie on the unit 
circle. For purely anti- ferromagnetic coupling, {za < 1, Zb < 1), they are on the imaginary 
axis whereas the "mixed" case turns out to be complicated and will not be discussed here. 

For simplicity, we will, from now on, stick to the purely ferromagnetic regime. If the chain 
were periodic, the zeros would be distributed homogeneously on the unit circle (they can be 
calculated analytically). In the thermodynamic limit, they fill the circle densely except near 
the point (1,0) on the real axis, where we are left with a gap since we do not have a phase 
transition at finite temperature (i.e., for < T < oo). One might perhaps expect the very 
same situation in the Fibonacci case, but the latter is always good for a surprise. 

In Fig. 1, we show the location of the zeros of Zn{za, Zb,w) in the complex w-plane 
for ferromagnetic couplings Za = 3/2, Zb = 100 (which corresponds to Ka ~ 0.203 and 
Kb ~ 2.306) and n = 8,9, 10. As expected, the zeros are clearly located on the unit circle, 
and there is still a large gap near (1,0) as it must be because we still cannot have a phase 
transition at finite temperature. However, an additional gap structure in the distribution of 
the zeros on the unit circle is apparent. It turns out that this gap structure does not depend 
on the actual values of the coupling constants (as long as they are still ferromagnetic), 
just the gap widths change and the gaps vanish if Za and Zb become equal which of course 
corresponds to the periodic case. The large difference between the two couplings used in 
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our pictures was chosen to show the effect clearly - for a small difference one might miss it, 
though it is still there! 



Figure 1: Zeros of the partition function Zn(za, Zb, w) of the Fibonacci Ising chain in the 
field variable w = exp{2h) for Za = 3/2, Zf, = 100, and n = 8, 9, 10. 
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In Fig. 2, we present the integrated density of the zeros on the unit circle (i.e., the 
integrated density of their angles (in units of 27r) with respect to the real axis) for the same 
set of parameter values as in Fig. 1. It turns out that this has exactly the structure which one 



would expect from the general gap labeling theorem of Bellissard and coworkers |I3| , 0, |T2 
which however applies to the integrated density of states (IDOS) of energy spectra of certain 
discrete Hamiltonians. There, the limit set of the IDOS values at the gaps has the form (see 
1^] and references therein) 

u V 
- + - 

T T'^ 



G z \ n [0, 1] (2.8) 



(/i,z/) : — > Z + Z^ (2-9) 



for the Fibonacci case, where r = is the golden ratio. The widest gaps in Fig. 2 have 
been labeled with the corresponding indices v) which for the finite system belong to the 
rational values 

— 1 n^oo ^ ^ 

fn+l T 

of the integrated density. Note that any two successive Fibonacci numbers and fn are 
coprime, hence every integer can be written as a linear combination of them with integral 
coefficients. Generically, it appears that at all the "allowed" values one actually observes 
gaps in the density distribution, i.e., "all gaps are open". The obvious symmetry in the 
pictures which relates the gaps with labels (/i, u) and (1 — /x, 1 — z/) stems from the reflection 
symmetry of the zero pattern with respect to the real axis which corresponds to a change of 
the sign in the magnetic field h. 
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Figure 2: Integrated density of tiie zeros in Fig. 1 on the unit circle. The corresponding 
gap labels (see Eq. (2.6)) for the widest gaps are also shown. 




Admittedly, we did not present a rigorous argument which explains why the gap labeling 
theorem for the IDOS of one-dimensional Schodinger operators [1^, |12] describes the 
distributions of our partition function zeros on the unit circle. However, the agreement is 
certainly convincing and suggests to think of the zeros (or rather of their arguments) as the 
eigenvalues of a Schrodinger type operator. This is also one relatively easy way to prove the 
Lee- Yang circle theorem (or the corresponding "line theorem" in the antiferromagnetic case) 
for the periodic ID Ising model. In this context, we also mention an old idea of Hilbert's, 
namely the possible connection between the imaginary parts of the (non-trivial) zeros of the 
Riemann (^-function and the eigenvalues of a - so far unknown - hermitian Hamiltonian. 
This has recently also been linked to some typical aspects of quantum chaos, see [15| and 
references therein. 



3 Ising model on the triangle tiling 

The ID case was presented for two main reasons. On the one hand, calculations are either 
possible analytically or can be made rigorous. On the other hand, a new phenomenon - the 
appearance of gaps - could be seen. 

Nevertheless, investigations of Ising models on graphs of higher dimension are in order. 
If one is not interested in approximative calculations (which we are not), one encounters true 
difficulties on the non-periodic ground. Even the existence of local deflation/inflation sym- 
metry does not seem to allow exact renormalization schemes for electronic models, compare 



[p^ , 0, and we were not able so far to find one for Ising type models, either. 

So, the best thing one can then do is the exact calculation of the partition sum for 
finite patches, followed by a numerical approach of the thermodynamic limit. Even this is a 
difficult task which requires a good choice of the graph, i.e., the non-periodic tiling. From 
our experience with Ising quantum chains in ID |]TB| - which can be seen as anisotropic 
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limits of 2D classical systems - we decided to take a quasiperiodic example in order to stay 
free of effects caused by uncontrollable fluctuations. 



Figure 3: Initial decagonal patch of the Tiibingen triangle pattern. 




The triangle pattern developed in Tiibingen |T^, ^ has the advantage that we can 
start from a patch with decagonal boundary that is - in a natural way - partitioned into 
sectors, each covering an angle of 27r/10, compare Fig. 3. This structure suggests the use of 
corner transfer matrices, compare |Q, to calculate the partition sum and the magnetization 
at the centre of the patch. However, the patch is not a repetition of ten equal sectors, 
wherefore we cannot reduce the problem to the CTM of one single sector. On the other 
hand, the CTM's of sectors of opposite orientation (indicated by the arrow in the inital 
patch) are the transposed matrices of each other. Giving those spins half weight which lie 
on the boundaries of the sectors, we can obtain the partition function as 



tr 



((M^M^M^) ■ (M^M^M^)*) . (3.1) 



This can then be evaluated again by algebraic manipulation packages. 

In what follows, we use the same notation as in the Fibonacci case above but with indices 
s (for short) and I (for long) in place of a and h. In particular, we also assume that the 
magnetic field is uniform. Here, Zn{zsi z^, w) now denotes the partition function of the Ising 
model on the patch which is obtained by n deflation steps (see Fig. 4) from the initial 
decagonal patch shown in Fig. 3 (which corresponds to n = 0). It is impossible to define 
periodic boundary conditions on our tiling without destroying the tenfold symmetry and the 
sector structure, wherefore we either use fixed boundary conditions (i.e., all spins on the 
decagonal boundary of the patch are frozen to be +1) or free boundary conditions (i.e., we 
sum over all values of the boundary spins). 
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Figure 4: Deflation rules for the Tiibingen triangle pattern. 




Figure 5: One sector of Fig. 3 after n deflation steps (see Fig. 4) for n < 4. The dashed 
line indicates a sector that also leads to a patch with decagonal boundary, 
which however is not obtained by entire deflation steps from the initial patch 
shown in Fig. 3. 




M=4 



For fixed boundary conditions, the zeros of tlie partition function Zn{zs, ze, w) (for finite 
n) in tlie field variable w are no longer located on the unit circle, as this case (in contrast 
to the free boundary case) is not covered by the circle theorem of Lee and Yang This 
is clearly seen in Fig. 6 which shows the locations of the zeros in the complex w-plane for 
couplings Zg = Zi = 4 and three different patch sizes. It is plausible that the zeros approach 
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the unit circle in the thermodynamic hmit. However, the angular distribution of the zeros 
is remarkably regular, there is no apparent gap structure as in the one- dimensional case. 
We also computed the zeros for different couplings Zg and for fixed and free boundary 
conditions which show the same behaviour. 



Figure 6: Zeros of the partition function Z„(4,4, w) (fixed boundary conditions) in the 
field variable w for n — 1,2,3. 
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Let us also have a look at the zeros in the other variables. Here, we restrict ourselves 
to the zero-field partition function Zn{zs,Zf^, 1) and use fixed boundary conditions. For the 
"isotropic" case z = Zg = z^, the zeros of ^3(2;,^, 1) in the variable z are shown in Fig. 7. 
Although the distribution of the zeros in the complex z plane is apparently not simple (i.e., 
the zeros do not appear to lie on simple curves), the zeros close to the real axis contain 
information about the location of the critical point and, in principle, also about the critical 
exponents ||^. In Table 1, we give the numerical values of the zero closest to the real line 
with a real part greater than one for the following five choices of (z^, Z(): (1, 2;), (2;, z^), (z, z), 
{z'^,z), and {z, 1), i.e., in addition to the "isotropic" case we look at those cases where one 
of the coupling constants vanishes and where one coupling constant is twice as large as the 
other. 
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Figure 7: Zeros of the partition function ^3(2, z, 1) (fixed boundary conditions) for the 
"isotropic" zero-field case (equal coupling for long and short bonds) in the 
temperature variable z. 
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Table 1: Partition function zero located closest to the real line and with a real part 
greater than one, for vanishing magnetic field and fixed boundary conditions. 



n 


2;^ = 1 


Zs — Zj^ 


Zs = Zl 


zl = Zl 


Z(. = \ 


f 


f.35f4±0.9560i 


1.2396 ±0.2819i 


1.3249 ±0.4315i 


1.1987± 0.2492i 


1.4511 ± 1.3404i 


2 


f.8056±0.8702i 


1.3270 ±0.2006i 


1.4709 ±0.3238i 


1.2788 ±0.1844i 


2.4447 ± 1.3460i 


3 


1.9089 ±0.6472i 


1.3670 ±0.1412i 


1.5249 ±0.2322i 


1.2983 ±0.1320i 


3.0508 ± 1.1471i 



The approximants for the critical couphngs are in good agreement with our values ob- 
tained from the behaviour of the specific heat and the centre spin magnetization of the model 
(with the exception of ze = 1 where the centre spin is isolated), but we omit details here. 
Let us only remark that the nature of the critical point looks very much like that of the 
periodic case though further calculations are needed to approve that. 
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4 Concluding Remarks 



The distribution of partition function zeros on the unit circle shows a gap structure for 
ferromagnetic Ising models on aperiodic chains. Although we have only demonstrated this 
phenomenon for the ubiquitious Fibonacci chain, one can translate each step to any of the 
other chains obtained by substitution rules - no matter whether one restricts oneself to the 
two-letter case or not. 

The 2D case did not show any apparent gap structure for the magnetic field zeros (and 
thus resembles the situation of the electronic spectra again). The partition function zeros in 
the temperature variable do not seem to lie on simple curves, even in the "isotropic" case 
where the couplings for short and long bonds are identical - more work is to be done to 
clarify this point. On the other hand, one can clearly see the existence of a phase transition 
because the zeros "pinch" the real axis. For given (finite) coupling constants, the critical 
point has finite T^.. 

Now, looking at the limiting cases Zg = 1 oi = 1, one gets the impression that heads 
to 0, quite similarly as in the case of the square lattice 0. However, the corresponding 
graph with one type of bonds "switched off" is by no means one-dimensional (or at least 
not in an obvious way) - wherefore the question arises what this means. To clear this up, 
one should treat first a model where the location of the critical point can be found exactly, 
without numerical estimates. This is indeed possible for another quasiperiodic pattern, the 
so-called Labyrinth tiling. The same phenomenom occurs there, and we will discuss details 
elsewhere . 
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